%   CLASS GPS_Time
% =========================================================================
%
% DESCRIPTION
%   Class to manage times and dates in various format (GPS / UTC/ ...)
%
% EXAMPLE
%   settings = GPS_Time();
%
% CONSTRUCTOR SYNTAX
%   t = GPS_Time(matlab_time, <[]>, <is_gps = 1>, <0>);
%   t = GPS_Time(uint32(unix_time), fraction_of_second, <is_gps = 1>, <1>);
%   t = GPS_Time(time_matlab_reference, time_difference, <is_gps = 1>, <2>);
%
% FOR A LIST OF CONSTANTs and METHODS use doc GPS_Time
%
% COMMENTS
% The class stores arrays of time, not just a single element,
% it has been designed this way because MATLAB works best on arrays
%
% This is probably an overcomplicated class, but storing times in multiple
% ways allows speed improvements.

%--- * --. --- --. .--. ... * ---------------------------------------------
%               ___ ___ ___
%     __ _ ___ / __| _ | __|
%    / _` / _ \ (_ |  _|__ \
%    \__, \___/\___|_| |___/
%    |___/                    v 0.5.2 beta 1
%
%--------------------------------------------------------------------------
%  Copyright (C) 2009-2017 Mirko Reguzzoni, Eugenio Realini
%  Written by:       Gatti Andrea
%  Contributors:     Gatti Andrea, ...
%  A list of all the historical goGPS contributors is in CREDITS.nfo
%--------------------------------------------------------------------------
%
%   This program is free software: you can redistribute it and/or modify
%   it under the terms of the GNU General Public License as published by
%   the Free Software Foundation, either version 3 of the License, or
%   (at your option) any later version.
%
%   This program is distributed in the hope that it will be useful,
%   but WITHOUT ANY WARRANTY; without even the implied warranty of
%   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%   GNU General Public License for more details.
%
%   You should have received a copy of the GNU General Public License
%   along with this program.  If not, see <http://www.gnu.org/licenses/>.
%
%--------------------------------------------------------------------------
% 01100111 01101111 01000111 01010000 01010011
%--------------------------------------------------------------------------


classdef GPS_Time < handle
    
    properties (Constant, GetAccess = public)
        DAYS_IN_WEEK = uint32(7);               % Number of days in a week
        SEC_IN_DAY  = uint32(86400);            % Number of seconds in a day
        SEC_IN_HALF_WEEK = uint32(302400);      % Number of seconds in a half a week
        SEC_IN_WEEK = uint32(604800);           % Number of seconds in a half a week
        UNIX_GPS_SEC_DIFF = uint32(315964800);  % Seconds of difference between UNIX time and GPS time
        
        DEFAULT_DATE_FORMAT = 'yyyy/mm/dd HH:MM:SS'; % String representing the format of visualization of the time
        
        % Epochs in (UTC time), matlab format when a leap second happened
        LEAP_DATES_UTC =   ([ 723728; 724093; 724458; 725189; 726103;
            726834; 727199; 727746; 728111; 728476;
            729025; 729572; 730121; 732678; 733774;
            735051; 736146; 736696]);
        
        % datenum ignores cycle slips -> in UTC leap seconds happen "before"
        LEAP_DATES_GPS =   ... %LEAP_DATES_UTC + (1:(length(LEAP_DATES_UTC)))'/86400;
            ([ 723728.000011574; 724093.000023148; 724458.000034722; 725189.000046296; 726103.00005787; ...
            726834.000069444; 727199.000081019; 727746.000092593; 728111.000104167; 728476.000115741; ...
            729025.000127315; 729572.000138889; 730121.000150463; 732678.000162037; 733774.000173611; ...
            735051.000185185; 736146.000196759; 736696.000208333; ]);
        
        %       Generated by:
        %                      datenum( '1981/07/01'; '1982/07/01'; '1983/07/01'; '1985/07/01'; '1988/01/01'; ...
        %                               '1990/01/01'; '1991/01/01'; '1992/07/01'; '1993/07/01'; '1994/07/01'; ...
        %                               '1996/01/01'; '1997/07/01'; '1999/01/01'; '2006/01/01'; '2009/01/01'; ...
        %                               '2012/07/01'; '2015/07/01'; '2017/01/01'})
        
        GPS_ZERO = 723186;              % datenum('Jan 6, 1980')
        UNIX_ZERO = 719529;             % datenum('Jan 1, 1970')
    end
    
    properties (Constant, GetAccess = private)
        MAT_TIME = 0;       % time_type value for times stored in matalb format: (double array) time in days since January 1, 0000 (matlab datenum fomat), precision up to the 0.1 milliseconds
        UNIX_TIME = 1;      % time_type value for times stored in unix format: (uint32 array) time in seconds since January 1, 1970 (UNIX standard) + (double array) fraction of seconds
        REF_TIME = 2;       % time_type value for times stored in time_ref: (double) origin of the "time system" expressed in datenum format + (double array) difference in seconds w.r.t. time_ref
    end
    
    properties (SetAccess = private, GetAccess = public)
        log = Logger.getInstance(); % Handler to the log object
        
        time_type = 0;      % flag depending on its value different representation of time are possible
        
        % time_type == 0 MATLAB_TIME DEFAULT it supports up to ~0.1 ms precision
        
        mat_time            % (double array) time in days since January 1, 0000 (matlab datenum fomat), precision up to the 0.1 milliseconds
        
        % time_type == 1 UNIX_TIME it supports ps precision
        
        unix_time           % (uint32 array) time in seconds since January 1, 1970 (UNIX standard)
        unix_time_f         % (double array) fraction of seconds
        
        % time_type == 2 REFERENCED_TIME it supports ps precision
        
        time_ref            % (double) origin of the "time system" expressed in datenum format
        time_diff           % (double array) difference in seconds w.r.t. time_ref
        
        is_gps = true;       % define whether is gps time or UTC
        
        
        date_format = GPS_Time.DEFAULT_DATE_FORMAT;
        leap_seconds = 999;
    end
    
    % Function to simulate polymorphism
    methods (Access = 'public')
        function this = GPS_Time_str(this, string_time, is_gps)
            % Private constructor - simulate polymorphism - GPS_Time_str(string_time, is_gps)
            if (nargin == 3)
                if isempty(is_gps)
                    is_gps = true;
                end
                this.is_gps = is_gps;
            end
            date = sscanf(string_time,'%f%f%f%f%f%f')';
            if (date(1) < 80), date(1) = date(1) + 2000; end
            this.GPS_Time_mat(datenum(date), is_gps);
        end
        
        function this = GPS_Time_mat(this, matlab_time, is_gps)
            % Private constructor - simulate polymorphism - GPS_Time_mat(matlab_time, is_gps)
            this.time_type = 0;
            this.mat_time = matlab_time;
            if (nargin == 3)
                if isempty(is_gps)
                    is_gps = true;
                end
                this.is_gps = is_gps;
            end
        end
        
        function this = GPS_Time_unix(this, unix_time, fraction_of_second, is_gps)
            % Private constructor - simulate polymorphism - GPS_Time_unix(uint32(unix_time), fraction_of_second, is_gps)
            this.time_type = 1;
            this.unix_time = uint32(unix_time);
            this.unix_time_f = fraction_of_second;
            if (nargin == 4)
                if isempty(is_gps)
                    is_gps = true;
                end
                this.is_gps = is_gps;
            end
        end
        
        function this = GPS_Time_6col(this, date, is_gps)
            % Private constructor - simulate polymorphism - GPS_Time_6col(data_6col, is_gps)
            
            % check date format 2/4 digits
            date(date(:,1) < 70,1) = date(date(:,1) < 70,1) + 2000;
            date(date(:,1) > 2070,1) = date(date(:,1) < 70,1) - 100;
            this.time_type = 1;
            if (nargin == 3)
                if isempty(is_gps)
                    is_gps = true;
                end
            elseif (nargin == 2)
                is_gps = true;
            end
            
            unx_ref = 719529; % this.UNIX_ZERO
            
            % number of days since the beginning of Unix time
            % deltat   = (datenum([date(:,1), date(:,2), date(:,3)]) - unx_ref);
            % hack: datenummmx is faster cause it does not check argins
            unix_time   = uint32((datenummx(date(:,1:3)) - unx_ref)) * 86400; %#ok<PROPLC>
            
            if (size(date, 2) == 6)
                % I want to keep precision (let's use Unix time)
                s = floor(date(:,6));
                unix_time = unix_time + uint32(date(:,4)) * 3600 + uint32(date(:,5)) * 60 + uint32(s);  %#ok<PROPLC>
                fraction_of_second = round((date(:,6) - s)*1e14)/1e14;
            else
                fraction_of_second = 0;
            end
            
            this.is_gps = is_gps;
            this.unix_time = unix_time; %#ok<PROPLC>
            this.unix_time_f = fraction_of_second;
        end
        
        function this = GPS_Time_Week_Dow(this, week, dow)
            % Private constructor - simulate polymorphism - GPS_Time_Week_Dow(week, dow)
            this.unix_time = uint32((this.GPS_ZERO - this.UNIX_ZERO)*86400 + week*86400*7 + dow*86400);
            this.unix_time_f = zeros(size(this.unix_time));
            this.time_type = 1;
        end
        
        function this = GPS_Time_ref(this, time_matlab_reference, time_difference, is_gps)
            % Private constructor - simulate polymorphism - GPS_Time_mat(time_matlab_reference, time_difference, is_gps)
            this.time_type = 2;
            this.time_ref = time_matlab_reference;
            this.time_diff = time_difference(:);
            if (nargin == 4)
                if isempty(is_gps)
                    is_gps = true;
                end
                this.is_gps = is_gps;
            end
        end
        
        function this = appendStrTime(this, string_time, is_gps)
            % Append elements - appendMatTime(matlab_time, is_gps)
            if (nargin == 3)
                if isempty(is_gps)
                    is_gps = this.is_gps;
                end
            elseif (nargin == 2)
                is_gps = this.is_gps;
            end
            date = sscanf(string_time,'%f%f%f%f%f%f')';
            if (date(1) < 80), date(1) = date(1) + 2000; end
            this.append(GPS_Time(datenum(date), [], is_gps, 0));
        end
        
        function this = appendMatTime(this, matlab_time, is_gps)
            % Append elements - appendMatTime(matlab_time, is_gps)
            if (nargin == 3)
                if isempty(is_gps)
                    is_gps = this.is_gps;
                end
            elseif (nargin == 2)
                is_gps = this.is_gps;
            end
            this.append(GPS_Time(matlab_time, [], is_gps, 0));
        end
        
        function this = appendUnixTime(this, unix_time, fraction_of_second, is_gps)
            % Append elements -  appendUnixTime(uint32(unix_time), fraction_of_second, is_gps)
            if (nargin == 4)
                if isempty(is_gps)
                    is_gps = this.is_gps;
                end
            elseif (nargin == 3)
                is_gps = this.is_gps;
            end
            this.append(GPS_Time(unix_time, fraction_of_second, is_gps, 1));
        end
        
        function this = append6ColDate(this, date, is_gps)
            % Append elements -  this.append6ColDate(date, is_gps)
            
            if (nargin == 3)
                if isempty(is_gps)
                    is_gps = this.is_gps;
                end
            elseif (nargin == 2)
                is_gps = this.is_gps;
            end
            
            unx_ref = 719529; % this.UNIX_ZERO
            
            % number of days since the beginning of Unix time
            % deltat   = (datenum([date(:,1), date(:,2), date(:,3)]) - unx_ref);
            % hack: datenummmx is faster cause it does not check argins
            unix_time   = uint32((datenummx(date(:,1:3)) - unx_ref) * 86400); %#ok<PROPLC>
            
            if (size(date, 2) == 6)
                % I want to keep precision (let's use Unix time)
                this.toUnixTime();
                s = floor(date(:,6));
                unix_time = unix_time + uint32(date(:,4)) * 3600 + uint32(date(:,5)) * 60 + uint32(s);  %#ok<PROPLC>
                fraction_of_second = round((date(:,6) - s)*1e14)/1e14;
            else
                fraction_of_second = 0;
            end
            
            this.append(GPS_Time(unix_time, fraction_of_second, is_gps, 1)); %#ok<PROPLC>
        end
        
        function this = appendRefTime(this, time_matlab_reference, time_difference, is_gps)
            % Append elements - appendRefTime(time_matlab_reference, time_difference, is_gps)
            if (nargin == 4)
                if isempty(is_gps)
                    is_gps = this.is_gps;
                end
            elseif (nargin == 3)
                is_gps = this.is_gps;
            end
            
            this.append(GPS_Time(time_matlab_reference, time_difference, is_gps, 2));
        end
    end
    
    methods
        function this = GPS_Time( arg1, arg2, arg3, arg4)
            % Constructor
            % SYNTAX:
            %   t = GPS_Time(matlab_time, <[]>, <is_gps = 1>, <0>);
            %   t = GPS_Time(uint32(unix_time), fraction_of_second, <is_gps = 1>, <1>);
            %   t = GPS_Time(time_matlab_reference, time_difference, <is_gps = 1>, <2>);
            
            switch nargin
                % With one parameter I assume to read GPS time in matlab format (string or number)
                % one value (in seconds) with milliseconds precision since 1 Jan 0000
                case 1
                    if isa(arg1,'char')
                        % string time (matlab time)
                        this.GPS_Time_str(arg1);
                    else
                        if size(arg1,2) == 1
                            % matlab time
                            this.GPS_Time_mat(arg1);
                        else
                            this.GPS_Time_6col(arg1);
                        end
                    end
                    % With two parameters it can be a unix o ref time format
                case 2
                    if isa(arg1,'uint32')
                        % UNIX time
                        this.GPS_Time_unix(arg1, arg2);
                    else
                        % Ref Time
                        this.GPS_Time_ref(arg1, arg2);
                    end
                    % With three parameters the thirsd is a flag (is_gps)
                case 3
                    if isempty(arg2)
                        if isa(arg1,'char')
                            % string time (matlab time)
                            this.GPS_Time_str(arg1, arg3);
                        else
                            if size(arg1,2) == 1
                                % matlab time
                                this.GPS_Time_mat(arg1, arg3);
                            else
                                this.GPS_Time_6col(arg1, arg3);
                            end
                        end
                    else
                        if isa(arg1,'uint32')
                            % UNIX time
                            this.GPS_Time_unix(arg1, arg2, arg3);
                        else
                            % Ref Time
                            this.GPS_Time_ref(arg1, arg2, arg3);
                        end
                    end
                case 4
                    switch arg4
                        case 0 % GPS_Time.MAT_TIME
                            this.GPS_Time_mat(arg1, arg3);
                        case 1 % GPS_Time.UNIX_TIME
                            this.GPS_Time_unix(arg1, arg2, arg3);
                        case 2 % GPS_Time.REF_TIME
                            this.GPS_Time_ref(arg1, arg2, arg3);
                        case 3 % GPS_Time_Week_Dow (to unix time)
                            this.GPS_Time_Week_Dow(arg1, arg2);
                        otherwise
                            this.log.addError('Unrecognized time format!!!');
                    end
            end
            %this.computeLeapSeconds();
        end
        
        function addEpoch(this, arg1, arg2, arg3, arg4)
            % Add epochs to the object as matlab/unix/ref time
            % SYNTAX:
            %   t = t.appendEpoch(matlab_time, <[]>, <is_gps = 1>, <0>);
            %   t = t.appendEpoch(uint32(unix_time), fraction_of_second, <is_gps = 1>, <1>);
            %   t = t.appendEpoch(time_matlab_reference, time_difference, <is_gps = 1>, <2>);
            
            switch nargin
                % With one parameter I assume to read GPS time in matlab format
                % one value (in seconds) with milliseconds precision since 1 Jan 0000
                case 2
                    if isa(arg1,'char')
                        % string time (matlab time)
                        this.appendStrTime(arg1);
                    else
                        % matlab time
                        this.appendMatTime(arg1);
                    end
                    
                    % With two parameters it can be a unix o ref time format
                case 3
                    if isa(uint32(arg1),'uint32')
                        % UNIX time
                        this.appendUnixTime(arg1, arg2);
                    else
                        % Ref Time
                        this.appendRefTime(arg1, arg2);
                    end
                    % With three parameters the thirsd is a flag (is_gps)
                case 4
                    if isempty(arg2)
                        if isa(arg1,'char')
                            % string time (matlab time)
                            this.appendStrTime(arg1, arg3);
                        else
                            % matlab time
                            this.appendMatTime(arg1, arg3);
                        end
                    else
                        if isa(uint32(arg1),'uint32')
                            % UNIX time
                            this.appendUnixTime(arg1, arg2, arg3);
                        else
                            % Ref Time
                            this.appendRefTime(arg1, arg2, arg3);
                        end
                    end
                case 5
                    switch arg4
                        case 0 % GPS_Time.MAT_TIME
                            this.appendMatTime(arg1, arg3);
                        case 1 % GPS_Time.UNIX_TIME
                            this.appendUnixTime(arg1, arg2, arg3);
                        case 2 % GPS_Time.REF_TIME
                            this.appendRefTime(arg1, arg2, arg3);
                        otherwise
                            this.log.addError('Unrecognized time format!!!');
                    end
            end
            %this.computeLeapSeconds();
        end
        
        function import(this, time)
            % Copy from an object of the same type (alias to copyFrom)
            this.copyFrom(time)
        end
        
        function copyFrom(this, time)
            % Copy from an object of the same type
            this.time_type = time.time_type;
            this.mat_time = time.mat_time;
            this.unix_time = time.unix_time;
            this.unix_time_f = time.unix_time_f;
            this.time_ref = time.time_ref;
            this.time_diff = time.time_diff;
            this.is_gps = time.is_gps;
            this.date_format = time.date_format;
            this.leap_seconds = time.leap_seconds;
        end
        
        function copy = getCopy(this)
            % Get a copy of this
            copy = GPS_Time();
            copy.copyFrom(this);
        end
        
        function this = append(this, time, time_type, is_gps)
            % Append a GPS_Time object into the this
            
            % check if the object is empty
            if isempty(this.time_type)
                this.import(time); % copy time in this
            else
                
                % When not specified the new format is the format of the append
                if (nargin < 3)
                    time_type = time.time_type;
                end
                if (nargin < 4)
                    is_gps = time.is_gps;
                end
                
                % Merge is done in GPS time -> it has no ambiguity
                this.toGps;
                time.toGps;
                
                switch time_type
                    case 0 % I'm in MAT TIME
                        this.toMatlabTime();
                        time.toMatlabTime();
                        
                        this.mat_time = [this.mat_time(:); time.mat_time(:)];
                        
                    case 1 % I'm in UNIX TIME
                        this.toUnixTime();
                        time.toUnixTime();
                        
                        this.unix_time = [this.unix_time; time.unix_time];
                        this.unix_time_f = [this.unix_time_f; time.unix_time_f];
                        
                    case 2 % I'm in REF TIME
                        this.toRefTime();
                        time.toRefTime();
                        
                        % keep the reference of the this
                        this.time_diff = [this.time_diff; (time.time_diff + (this.time_ref - time.time_ref) * 86400)];
                end
                
                % Convert into correct format time
                if is_gps
                    this.toGps();
                else
                    this.toUtc();
                end
                
                % Compute leap seconds, if the original data stored them
                if (time.leap_seconds < 999) ||  (this.leap_seconds < 999)
                    this.computeLeapSeconds();
                end
            end
        end
    end
    
    % =========================================================================
    %    CONVERSIONS
    % =========================================================================
    
    methods (Access = 'private')
        function leap_seconds = computeLeapSeconds(this)
            % compute the number of leap seconds to subtract to GPS Time to obtain the UTC time
            if this.is_gps
                LEAP_DATES = this.LEAP_DATES_GPS;
            else
                LEAP_DATES = this.LEAP_DATES_UTC;
            end
            
            switch this.time_type
                case 0 % I'm in MAT TIME
                    if (numel(this.mat_time) == 1)
                        leap_seconds = find(this.mat_time > LEAP_DATES, 1, 'last');
                    else
                        leap_seconds = [find(this.mat_time(1) > LEAP_DATES, 1, 'last'); ...
                            find(this.mat_time(end) > LEAP_DATES, 1, 'last') ];
                        if (diff(leap_seconds) == 0) % the same leap_second for the whole period
                            leap_seconds = leap_seconds(1);
                        elseif (diff(leap_seconds) == 1) % there's one leap in the dataset
                            % search for the latest value before the second leap second
                            % max precision of leap dates is 1e-9
                            id_leap = find(round(this.mat_time*1e9) >= round(LEAP_DATES(leap_seconds(2))*1e9),1, 'first') - 1;
                            leap_seconds = [ones(id_leap,1) * leap_seconds(1); ones(numel(this.mat_time) - id_leap,1) * leap_seconds(2)];
                        else % there are multiple leaps in the dataset
                            lp = length(length(this.mat_time));
                            lp(1) = leap_seconds(1);
                            lp(end) = leap_seconds(2);
                            for i = 2:length(this.mat_time)-1
                                lp(i) = find(this.mat_time > LEAP_DATES, 1, 'last');
                            end
                            leap_seconds = lp;
                        end
                    end
                case 1 % I'm in UNIX TIME
                    if (numel(this.mat_time) == 1)
                        leap_seconds = find((double(this.unix_time)) / 86400 > (LEAP_DATES - 719529), 1, 'last');
                    else
                        leap_seconds = [ find((double(this.unix_time(1))) / 86400 > (LEAP_DATES - 719529), 1, 'last'); ...
                            find((double(this.unix_time(end))) / 86400 > (LEAP_DATES - 719529), 1, 'last')];
                        if (diff(leap_seconds) == 0) % the same leap_second for the whole period
                            leap_seconds = leap_seconds(1);
                        elseif (diff(leap_seconds) == 1)  % there's one leap in the dataset
                            % search for the latest value before the second leap second
                            id_leap = find((this.unix_time) >= uint32((LEAP_DATES(leap_seconds(2)) - 719529) * 86400),1, 'first') - 1;
                            leap_seconds = [ones(id_leap,1) * leap_seconds(1); ones(numel(this.unix_time) - id_leap,1) * leap_seconds(2)];
                        else % there are multiple leaps in the dataset
                            lp = length(length(this.unix_time));
                            lp(1) = leap_seconds(1);
                            lp(end) = leap_seconds(2);
                            for i = 2:length(this.unix_time)-1
                                lp(i) = find((double(this.unix_time(i))) / 86400 > (LEAP_DATES - 719529), 1, 'last');
                            end
                            leap_seconds = lp;
                        end
                    end
                    
                case 2 % I'm in REF TIME
                    if (numel(this.mat_time) == 1)
                        leap_seconds = find(this.time_ref + this.time_diff  > LEAP_DATES, 1, 'last');
                    else
                        leap_seconds = [find((this.time_ref + this.time_diff(1)) > LEAP_DATES, 1, 'last'); ...
                            find((this.time_ref + this.time_diff(end)) > LEAP_DATES, 1, 'last') ];
                        if (diff(leap_seconds) == 0) % the same leap_second for the whole period
                            leap_seconds = leap_seconds(1);
                        elseif (diff(leap_seconds) == 1) % there's one leap in the dataset
                            % search for the latest value before the second leap second
                            % max precision of leap dates is 1e-9
                            id_leap = find(round(this.time_diff*1e9) >= round((LEAP_DATES(leap_seconds(2)) - this.time_ref)*1e9),1, 'first') - 1;
                            leap_seconds = [ones(id_leap,1) * leap_seconds(1); ones(numel(this.time_diff) - id_leap,1) * leap_seconds(2)];
                        else % there are multiple leaps in the dataset
                            lp = length(length(this.time_ref));
                            lp(1) = leap_seconds(1);
                            lp(end) = leap_seconds(2);
                            for i = 2:length(this.time_ref)-1
                                lp(i) = find(this.time_ref + this.time_diff  > LEAP_DATES, 1, 'last');
                            end
                            leap_seconds = lp;
                        end
                    end
                    
            end
            if ~(this.is_gps)
                leap_seconds = -leap_seconds;
            end
            this.leap_seconds = leap_seconds;
        end
    end
    
    methods
        function toMatlabTime(this)
            % Convert the internal structure to Matlab Time, precision up to the 0.1 milliseconds precision
            switch this.time_type
                case 0 % I'm already in MAT TIME
                    % do nothing
                case 1 % I'm in UNIX TIME
                    this.time_type = 0;
                    % constants in matlab are slower than copied values :-( switching to values
                    % this.mat_time = double(this.unix_time) / this.SEC_IN_DAY + this.UNIX_ZERO + this.unix_time_f;
                    this.mat_time = ((double(this.unix_time) + this.unix_time_f) / 86400 + 719529);
                    this.unix_time = [];
                    this.unix_time_f = [];
                case 2 % I'm in REF TIME
                    this.time_type = 0;
                    this.mat_time = this.time_ref + this.time_diff / 86400;
                    this.time_ref = [];
                    this.time_diff = [];
            end
        end
        
        function toUnixTime(this)
            % Convert the internal structure to Unix Time, precision up to the ps precision
            switch this.time_type
                case 0 % I'm in MAT TIME
                    this.time_type = 1;
                    % constants in matlab are slower than copied values :-( switching to values
                    % time_s = (this.mat_time - this.UNIX_ZERO) * this.SEC_IN_DAY; % convert mat_time in seconds
                    % due to numerical error propagation I can keep only 4 decimal
                    time_s =  (this.mat_time - 719529) * 86400; % convert mat_time in seconds
                    this.unix_time = uint32(fix(round(time_s * 1e4) / 1e4));
                    rounding = 10.^(round(16 - max(0,log10(this.mat_time * 86400 + eps(this.mat_time * 86400))))); % 52 bits of mantissa
                    this.unix_time_f = round((time_s - double(this.unix_time)) .* rounding) ./ rounding;                    
                    second_correction = - floor(this.unix_time_f);
                    this.unix_time = this.unix_time - uint32(second_correction);
                    this.unix_time_f = this.unix_time_f + second_correction;
                    clear time_s
                    this.mat_time = [];
                case 1 % I'm already in UNIX TIME
                    % do nothing
                case 2 % I'm in REF TIME
                    this.time_type = 1;
                    % constants in matlab are slower than copied values :-( switching to values
                    % time_s = (this.mat_time - this.UNIX_ZERO) * this.SEC_IN_DAY; % convert mat_time in seconds
                    time_s = (this.time_ref - 719529) * 86400;
                    this.unix_time = uint32(fix(round((time_s + this.time_diff) * 1e4) / 1e4));
                    rounding = 2.^(round(47 - max(0,log2(time_s + eps(time_s))))); % 52 bits of mantissa -> but I need to consider only 47 bits to keep precision
                    this.unix_time_f = round((time_s - double(this.unix_time)) .* rounding) ./ rounding + this.time_diff;
                    second_correction = - floor(this.unix_time_f);
                    this.unix_time = this.unix_time - uint32(second_correction);
                    this.unix_time_f = this.unix_time_f + second_correction;
                    clear time_s
                    this.time_ref = [];
                    this.time_diff = [];
            end
        end
        
        function toRefTime(this)
            % Convert the internal structure to Reference Time, precision up to the ps precision
            switch this.time_type
                case 0 % I'm in MAT TIME
                    this.time_type = 2;
                    this.time_ref = fix(this.mat_time(1));
                    % due to numerical error propagation I can keep only 4 decimal digits
                    this.time_diff = (this.mat_time - this.time_ref) * 86400;
                    this.mat_time = [];
                case 1 % I'm in UNIX TIME
                    this.time_type = 2;
                    % constants in matlab are slower than copied values :-( switching to values
                    % time_d = double(this.unix_time) / this.SEC_IN_DAY + this.UNIX_ZERO;
                    time_d = double(this.unix_time) / 86400 + 719529;
                    this.time_ref = fix(time_d(1));
                    this.time_diff = round((time_d - this.time_ref) * 86400) + this.unix_time_f;
                    this.unix_time = [];
                    this.unix_time_f = [];
                case 2 % I'm already in REF TIME
                    % do nothing
            end
        end
        
        function toUtc(this)
            % Transform the internal allocation in UTF format (corrects for cycle-slip
            if (this.is_gps == true)
                if (this.leap_seconds >= 999)
                    this.leap_seconds = this.computeLeapSeconds();
                end
                this.addSeconds(-this.leap_seconds);
                this.is_gps = false;
                this.leap_seconds = -this.leap_seconds;
            end
        end
        
        function toGps(this)
            % Transform the internal allocation in GPS format (corrects for cycle-slips)
            if (this.is_gps == false)
                if (this.leap_seconds >= 999)
                    this.leap_seconds = this.computeLeapSeconds();
                end
                this.addSeconds(-this.leap_seconds);
                this.is_gps = true;
                this.leap_seconds = -this.leap_seconds;
            end
        end
        
        function leap_seconds = getLeapSeconds(this)
            % get the number of leap seconds to subtract to GPS Time to obtain the UTC time
            if (this.leap_seconds >= 999)
                this.leap_seconds = this.computeLeapSeconds();
            end
            leap_seconds = this.leap_seconds();
        end
    end
    
    % =========================================================================
    %    GETTERS
    % =========================================================================
    
    methods
        function [is_utc] = isUTC(this)
            % Get the memorization type
            is_utc = ~this.is_gps;
        end
        
        function [is_gps] = isGPS(this)
            % Get the memorization type
            is_gps = this.is_gps;
        end
        
        function [len] = numel(this)
            % get number of epochs
            if isempty(this.time_type)
                len = 0;
            else
                switch this.time_type
                    case 0 % I'm in MAT TIME
                        len = numel(this.mat_time);
                    case 1 % I'm in UNIX TIME
                        len = numel(this.unix_time);
                    case 2 % I'm in REF TIME
                        len = numel(this.time_diff);
                end
            end
        end
        
        function n_element = length(this)
            % Get the total number of element stored in the object
            if isempty(this.time_type)
                n_element = 0;
            else
                switch this.time_type
                    case 0 % I'm in MAT TIME
                        n_element = length(this.mat_time);
                    case 1 % I'm in UNIX TIME
                        n_element = length(this.unix_time);
                    case 2 % I'm in REF TIME
                        n_element = length(this.time_diff);
                end
            end
        end
        
        function n_element = getLen(this)
            % Get the total number of element stored in the object
            n_element = this.length();
        end
        
        function [n_epochs] = getExpectedLen(this)
            % get the number of epochs (with constant rate) that the object should contain
            n_epochs = round((this.last.getMatlabTime() - this.first.getMatlabTime()) * (86400 / this.getRate())) + 1;
        end
        
        function [empty]  = isempty(this)
            % return the status of emptyness of the object
            switch this.time_type
                case 0 % I'm in MAT TIME
                    empty = isempty(this.mat_time);
                case 1 % I'm in UNIX TIME
                    empty = isempty(this.unix_time);
                case 2 % I'm in REF TIME
                    empty = isempty(this.time_diff);
            end
        end
        
        function [nan_list] = isnan(this)
            % return the logical NaN of the object
            switch this.time_type
                case 0 % I'm in MAT TIME
                    nan_list = isnan(this.mat_time);
                case 1 % I'm in UNIX TIME
                    nan_list = isnan(this.unix_time);
                case 2 % I'm in REF TIME
                    nan_list = isnan(this.time_diff);
            end
        end
        
        function [rate]  = getRate(this)
            % get observation rate approximated at 3 digits
            switch this.time_type
                case 0 % I'm already in MAT TIME
                    rate = round(median(diff(this.mat_time*86400)) * 1e3) / 1e3;
                case 1 % I'm in UNIX TIME
                    % constants in matlab are slower than copied values :-( switching to values
                    % this.mat_time = double(this.unix_time) / this.SEC_IN_DAY + this.UNIX_ZERO + this.unix_time_f;
                    tmp_time = double(this.unix_time) + this.unix_time_f;
                    rate = round(median(diff(tmp_time)) * 1e3) / 1e3;
                case 2 % I'm in REF TIME
                    rate = round(median(diff(this.time_diff / 86400)) * 1e3) / 1e3;
            end
        end
        
        function [mat_time] = getMatlabTime(this)
            % get Matlab Time, precision up to the 0.1 milliseconds precision
            switch this.time_type
                case 0 % I'm already in MAT TIME
                    mat_time = this.mat_time;
                case 1 % I'm in UNIX TIME
                    % constants in matlab are slower than copied values :-( switching to values
                    % this.mat_time = double(this.unix_time) / this.SEC_IN_DAY + this.UNIX_ZERO + this.unix_time_f;
                    mat_time = ((double(this.unix_time) + this.unix_time_f) / 86400 + 719529);
                case 2 % I'm in REF TIME
                    mat_time = this.time_ref + this.time_diff / 86400;
            end
        end
        
        function [unix_time, unix_time_f, second_correction] = getUnixTime(this)
            % get Unix Time, precision up to the ps precision
            second_correction = 0;
            switch this.time_type
                case 0 % I'm in MAT TIME
                    % constants in matlab are slower than copied values :-( switching to values
                    % time_s = (this.mat_time - this.UNIX_ZERO) * this.SEC_IN_DAY; % convert mat_time in seconds
                    % due to numerical error propagation I can keep only 4 decimal
                    time_s =  (this.mat_time - 719529) * 86400; % convert mat_time in seconds
                    unix_time = uint32(fix(round(time_s * 1e4) / 1e4));
                    unix_time_f = round((time_s - double(unix_time)) * 1e4) / 1e4;
                    second_correction = - floor(unix_time_f);
                    unix_time = unix_time - uint32(second_correction);
                    unix_time_f = unix_time_f + second_correction;
                case 1 % I'm already in UNIX TIME
                    unix_time = this.unix_time;
                    unix_time_f = this.unix_time_f;
                case 2 % I'm in REF TIME
                    % constants in matlab are slower than copied values :-( switching to values
                    % time_s = (this.mat_time - this.UNIX_ZERO) * this.SEC_IN_DAY; % convert mat_time in seconds
                    time_s = (this.time_ref - 719529) * 86400;
                    unix_time = uint32(fix(round((time_s + this.time_diff) * 1e4) / 1e4));
                    rounding = 2.^(round(47 - max(0,log2(time_s + eps(time_s))))); % 52 bits of mantissa -> but I need to consider only 47 bits to keep precision
                    unix_time_f = round((time_s - double(unix_time)) .* rounding) ./ rounding + this.time_diff;
                    second_correction = - floor(unix_time_f);
                    unix_time = unix_time - uint32(second_correction);
                    unix_time_f = unix_time_f + second_correction;
            end
        end
        
        function date6col = get6ColDate(this)
            % get time, 6 column format
            %
            % SYNTAX
            %   date6col = get6ColDate(this)
            time_base = this.getCopy; 
            time_base.toUnixTime; % to be safe
            time_base.unix_time_f = zeros(size(time_base.unix_time_f)); % removing fractional part
            date6col = datevec(time_base.getMatlabTime);
            date6col(:,6) = floor(date6col(:,6));
            tmp = GPS_Time(date6col);
            tmp.toUnixTime; % to be safe
            fractional_part = round((tmp - this) * 1e10) * 1e-10;
            date6col(:,6) = date6col(:,6) - fractional_part ;
        end
        
        function changeRef(this, new_time_mat_ref)
            % change the reference time in use by the object to store the data
            
            this.toRefTime(); % the object should store data into ref_time data type
            this.time_diff = this.time_diff + (this.time_ref - new_time_mat_ref) * 86400;
            this.time_ref = new_time_mat_ref;
        end
        
        function [time_diff, time_ref] = getRefTime(this, new_time_mat_ref)
            % get Reference Time, precision up to the ps precision
            % SYNTAX: [time_diff, time_ref] = this.getRefTime(new_time_mat_ref)
            switch this.time_type
                case 0 % I'm in MAT TIME
                    time_ref = fix(this.mat_time(1));
                    % due to numerical error propagation I can keep only 4 decimal digits
                    time_diff = (this.mat_time - time_ref) * 86400;
                    if (nargin == 2)
                        % optional change of reference
                        time_diff = this.time_diff + (this.time_ref - new_time_mat_ref) * 86400;
                        time_ref = new_time_mat_ref;
                    end
                case 1 % I'm in UNIX TIME
                    % constants in matlab are slower than copied values :-( switching to values                    
                    % this.UNIX_ZERO == 719529; 
                    if nargin == 2
                        time_ref = fix((new_time_mat_ref(1)  - 719529.0) * 86400);
                        rounding = 10.^(round(16 - max(0,log10(ceil(new_time_mat_ref(1)) * 86400 + eps(new_time_mat_ref(1) * 86400))))); % 52 bits of mantissa
                        time_ref_f = round(rem((new_time_mat_ref(1)  - 719529.0) * 86400, 1) * rounding) / rounding;
                    else
                        time_ref = this.unix_time(1);
                        time_ref_f = 0;
                    end
                    time_diff = double(int64(this.unix_time) - int64(time_ref)) + this.unix_time_f - time_ref_f;
                    time_ref = (time_ref / 86400 + 719529.0) + time_ref_f / 86400;
                case 2 % I'm in REF TIME
                    if (nargin == 2)
                        % optional change of reference
                        time_diff = this.time_diff + (this.time_ref - new_time_mat_ref) * 86400;
                        time_ref = new_time_mat_ref;
                    else
                        time_ref = this.time_ref;
                        time_diff = this.time_diff;
                    end
            end
            
        end
        
        function [gps_week, gps_sow, gps_dow] = getGpsWeek(this,id)
            % get Reference Time, precision up to the ps precision
            % SYNTAX: [gps_week, gps_sow, gps_dow] = this. getGpsWeek();
            gps_time = this.getCopy();
            gps_time.toGps();
            [unix_time, unix_time_f] = gps_time.getUnixTime(); %#ok<PROP>
            if nargin == 2
                unix_time = unix_time(id);
                unix_time_f = unix_time_f(id);
            end
            [gps_week, gps_sow, gps_dow] = gps_time.unixTimeToGps(unix_time, unix_time_f); %#ok<PROP>
        end
        
        function [gps_time] = getGpsTime(this, gps_offset)
            % Get time as number of seconds from Jan 6, 1980
            % SYNTAX: [gps_time] = this.getGpsTime(gps_offset)
            [unix_time, unix_time_f] = this.getUnixTime(); %#ok<PROPLC>
            if nargin == 1
                % gps_time = double(unix_time -  GPS_Time.UNIX_GPS_SEC_DIFF) + unix_time_f;
                gps_time = double(unix_time - uint32(315964800)) + unix_time_f;  %#ok<PROPLC>
            else
                max_digits = 15-ceil(log10(gps_offset));
                gps_time = round((double(unix_time - uint32(315964800)) - gps_offset) * max_digits) / max_digits + unix_time_f; %#ok<PROPLC>
            end
        end
        
        function [year, doy] = getDOY(this)
            % get Reference Time, precision up to the ps precision
            % SYNTAX: [year, doy] = getDOY(this)
            utc_time = this.getCopy();
            utc_time.toGps();
            utc_time.toMatlabTime();
            
            [year, ~] = datevec(utc_time.mat_time);
            doy = floor(utc_time.mat_time - datenum(year,1,1)) + 1; % days from the beginning of the year
        end
        
        function [year, month, day, hour, minute, second] = getCalEpoch(this,idx)
            % get year month doy hour minute second
            % SYNTAX: [year, month, doy, hour, minute, second] = getCalEpoch(this,idx)
            if nargin == 1
                idx = ones(this.length(),1) > 0;
            end
            str_time=this.toString();
            year = str2num(str_time(idx,1:4));
            month = str2num(str_time(idx,6:7));
            day = str2num(str_time(idx,9:10));
            hour = str2num(str_time(idx,12:13));
            minute = str2num(str_time(idx,15:16));
            second = str2num(str_time(idx,18:27));
        end
        
        function date_string = toString(this, date_format)
            % Convert a date to string format
            if this.isempty()
                date_string = '';
            else
                if (nargin == 2)
                    time = this.getMatlabTime();
                    time(isnan(time)) = 0;
                    date_string = datestr(time, date_format);
                else
                    time = this.getMatlabTime();
                    time(isnan(time)) = 0;
                    [~, fraction_of_seconds] = this.getUnixTime();
                    date_6col = datevec(time - (round(fraction_of_seconds)) / 86400);
                    date_6col(:,6) = round(date_6col(:,6)) + fraction_of_seconds;
                    date_string = reshape(sprintf('%04d/%02d/%02d %02d:%02d:%010.7f', date_6col')',27,numel(time))';
                end
                if (nargin == 1)
                    if this.isGPS()
                        date_string = [char(date_string(:,:)) char(repmat(' GPS',size(date_string,1),1))];
                    else
                        date_string = [char(date_string(:,:)) char(repmat(' UTC',size(date_string,1),1))];
                    end
                elseif strfind(date_format,'TTT') %#ok<*STRIFCND>
                    if this.isGPS()
                        date_string = reshape(regexprep(serialize(date_string')', 'TTT','GPS'), size(date_string,2), size(date_string,1))';
                    else
                        date_string = reshape(regexprep(serialize(date_string')', 'TTT','UTC'), size(date_string,2), size(date_string,1))';
                    end
                end
            end
        end
        
        function date_string = toStringGpsWeek(this)
            % Convert a date to string format
            if this.isempty()
                date_string = '';
            else
                [w, s, d] = this.getGpsWeek();
                date_string = sprintf('Week %d - dow %d - sow %d\n', w, d, s);
            end
        end
        
        function round(this, rounding)
            % round the time
            % SYNTAX: 
            %   this.round(rounding);
            % EXAMPLE:
            %   this.round(1e-5);
            switch this.time_type
                case 0 % I'm in MAT TIME
                    this.mat_time = (round(this.mat_time * 86400 ./ rounding) .* rounding) / 86400;
                case 1 % I'm in UNIX TIME
                    this.unix_time_f = round(this.unix_time_f ./ rounding) .* rounding;
                    idx =  this.unix_time_f >= 1;                    
                    this.unix_time(idx) = this.unix_time(idx) + uint32(ceil(this.unix_time_f(idx)));
                    this.unix_time_f(idx) = this.unix_time_f(idx) - ceil(this.unix_time_f(idx));
                    
                    idx =  this.unix_time_f < 0;                    
                    this.unix_time(idx) = this.unix_time(idx) - uint32(-floor(this.unix_time_f(idx)));
                    this.unix_time_f(idx) = this.unix_time_f(idx) - floor(this.unix_time_f(idx));
                case 2 % I'm in REF TIME
                    this.time_diff = round(this.time_diff ./ rounding) .* rounding;
            end
        end
        
        function new_obj = getEpoch(this, id)
            % Overloading of the operator index ()
            % get a copy of the obj containing only the selected epoch id of time
            % SYNTAX this.getEpoch(id)
            
            if islogical(id)
                max_id = find(id == true, 1, 'last');
            else
                max_id = max(id);
            end
            
            switch this.time_type
                case 0 % I'm in MAT TIME
                    if (length(this.mat_time) >= max_id)
                        new_obj = GPS_Time(this.mat_time(id), [], this.is_gps);
                    else
                        new_obj = GPS_Time();
                    end
                case 1 % I'm in UNIX TIME
                    if (length(this.unix_time) >= max_id)
                        new_obj = GPS_Time(uint32(this.unix_time(id)), this.unix_time_f(id), this.is_gps);
                    else
                        new_obj = GPS_Time();
                    end
                case 2 % I'm in REF TIME
                    if (length(this.time_diff) >= max_id)
                        new_obj = GPS_Time(this.time_ref, this.time_diff(id), this.is_gps);
                    else
                        new_obj = GPS_Time();
                    end
            end
        end
        function remEpoch(this, id)
            %DESCRIPTION: remove epochs by index
                    % SYNTAX this.remEpoch(id)
            
            if islogical(id)
                max_id = find(id == true, 1, 'last');
            else
                max_id = max(id);
            end
            
            switch this.time_type
                case 0 % I'm in MAT TIME
                    if (length(this.mat_time) >= max_id)
                        this.mat_time(id) = [];
                    end
                case 1 % I'm in UNIX TIME
                    if (length(this.unix_time) >= max_id)
                        this.unix_time(id) = [];
                        this.unix_time_f(id) = [];
                    end
                case 2 % I'm in REF TIME
                    if (length(this.time_diff) >= max_id)
                        this.time_diff(id) = [];
                    end
            end
        end

        function new_obj = first(this, id_subset)
            % Get first element stored in GPS_Time
            if (nargin == 1)
                new_obj = this.getEpoch(1);
            else
                if islogical(id_subset)
                    new_obj = this.getEpoch(find(id_subset, 1, 'first'));
                else
                    new_obj = this.getEpoch(id_subset(1));
                end
            end
        end
        
        function new_obj = last(this, id_subset)
            % Get last element stored in GPS_Time
            if (nargin == 1)
                new_obj = this.getEpoch(this.length());
            else
                if islogical(id_subset)
                    new_obj = this.getEpoch(find(id_subset(1 : this.length()), 1, 'last'));
                else
                    id_subset = id_subset(id_subset < this.length());
                    new_obj = this.getEpoch(id_subset(end));
                end
            end
        end
        function gps_time_subset = getSubSet(this, index)
            % create a copy of the object having only a subset of time_diff
            gps_time_subset = this.getCopy();
            switch  this.time_type
                case 0
                    gps_time_subset.mat_time = gps_time_subset.mat_time(index);
                case 1
                    gps_time_subset.unix_time = gps_time_subset.unix_time(index);
                    gps_time_subset.unix_time_f = gps_time_subset.unix_time_f(index);
                case 2
                    gps_time_subset.unix_time = gps_time_subset.unix_time(index);
            end
        end
    end
    
    % =========================================================================
    %    SETTERS
    % =========================================================================
    
    methods
        function setGPS(this, is_gps)
            % Set the memorization type
            this.is_gps = is_gps;
        end

        function setDateFormat(this, date_format)
            % Change the default value for "date format"
            this.date_format = date_format;
        end
        
        function delId(this, id)
            % Delete time with id = id
            switch this.time_type
                case 0 % I'm in MAT TIME
                    this.mat_time(id) =  [];
                case 1 % I'm in UNIX TIME
                    this.unix_time(id) = [];
                    this.unix_time_f(id) = [];
                case 2 % I'm in REF TIME
                    this.time_diff(id) = [];
            end
        end
        
        function delLast(this)
            % Delete last element stored in GPS_Time
            this.delId(this.length());
        end
    end
    
    % =========================================================================
    %    OPERATIONS
    % =========================================================================
    
    methods (Access = 'public')
        function addIntSeconds(this, n_seconds)
            % Add an integer number of seconds to all the times
            switch this.time_type
                case 0 % I'm in MAT TIME
                    this.mat_time = this.mat_time + n_seconds / 86400;
                case 1 % I'm in UNIX TIME
                    this.unix_time = uint32(int64(this.unix_time) + int64(n_seconds));
                case 2 % I'm in REF TIME
                    this.time_diff = this.time_diff + n_seconds;
            end
        end
        
        function addSeconds(this, n_seconds)
            % Add a floating point number of seconds to all the times
            switch this.time_type
                case 0 % I'm in MAT TIME
                    this.mat_time = this.mat_time + n_seconds / 86400;
                case 1 % I'm in UNIX TIME
                    this.unix_time = uint32(int64(this.unix_time)  + int64(fix(n_seconds)));
                    this.unix_time_f = this.unix_time_f + rem(n_seconds,1);
                    this.unix_time = uint32(int64(this.unix_time) + int64(floor(this.unix_time_f)));
                    this.unix_time_f = this.unix_time_f - floor(this.unix_time_f);
                case 2 % I'm in REF TIME
                    this.time_diff = this.time_diff + n_seconds;
            end
        end
        
        function res = lt(gt_1, gt_2)
            %%% DESCRIPTION: overload of '<' function
            %%% get unix time
            [unix_time1, unix_time_f1] = gt_1.getUnixTime();
            [unix_time2, unix_time_f2] = gt_2.getUnixTime();
            
            res = (unix_time1 < unix_time2) | ((unix_time1 == unix_time2)  & (unix_time_f1 < unix_time_f2)) ;
        end
        
        function res = gt(gt_1, gt_2)
            %%% DESCRIPTION: overload of '>' function
            %%% get unix time
            [unix_time1, unix_time_f1] = gt_1.getUnixTime();
            [unix_time2, unix_time_f2] = gt_2.getUnixTime();
            
            res = (unix_time1 > unix_time2) |( (unix_time1 == unix_time2) & (unix_time_f1 > unix_time_f2) );
        end
        
        function res = eq(gt_1, gt_2)
            %%% DESCRIPTION: check if two time are equals up to precision
            
            [~, sec, sec_f] = gt_1 -gt_2;
            prec = max(gt_1.getPrecision, gt_2.getPrecision);
            res = abs(sec+sec_f) < prec;
            
        end
        
        function [sec, sec_i, sec_f] = minus(gt_1, gt_2)
            % OUTPUT:
            %    sec = difference in seconds
            %    sec = unit part
            %    sec_f = fractional difference in seconds
            % DESCRIPTION: overload of '-' function
            %%% get unix time
            [unix_time1, unix_time_f1] = gt_1.getUnixTime();
            [unix_time2, unix_time_f2] = gt_2.getUnixTime();
            sec_i = int64(unix_time1) - int64(unix_time2);
            sec_f = unix_time_f1 - unix_time_f2;
            
            
            % make the two values consistent
            idx_sec = sec_i > 0;
            idx_sec_f = sec_f >0;
            idx_conflict = xor(idx_sec,idx_sec_f);
            if sum(idx_conflict) > 0
                sec_i(idx_conflict) = sec_i(idx_conflict) + int64(sign(sec_f(idx_conflict)));
                sec_f(idx_conflict) = sec_f(idx_conflict) - sign(sec_f(idx_conflict));
            end
            sec = double(sec_i) + sec_f;
            
        end
        
        function prec = getPrecision(this)
            % get precison of current time
            switch this.time_type
                case 0 % I'm in MAT TIME
                    prec = eps(this.mat_time);
                case 1 % I'm in UNIX TIME
                    prec = eps(this.unix_time_f);
                case 2 % I'm in REF TIME
                    prec = eps(this.time_diff);
            end
        end
    end
    
    % =========================================================================
    %    STATIC UNIX TIME
    % =========================================================================
    
    methods (Static, Access = 'public')
        function [unix_time, unix_time_f] = matToUnixTime(mat_time)
            % Conversion From MATLAB (expressed in days) type to UNIX TIME (expressed in uint32 seconds)
            time_s = round((mat_time - 719529) * 86400 * 1e4)/1e4; % convert mat_time in seconds
            
            if (nargout == 2)
                unix_time_f = rem(time_s,1);
            end
            unix_time = uint32(fix(time_s));
        end
        
        function [unix_time, unix_time_f] = dateStringToUnixTime(time_string)
            % Conversion From STRING type to UNIX TIME (expressed in uint32 seconds) max precision ~0.1 ms
            [unix_time, unix_time_f] = GPS_Time.matToUnixTime(datenum(time_string));
        end
        
        function [unix_time, unix_time_f] = gpsToUnixTime(gps_week, gps_time)
            % Conversion (shift) from GPS time (January 6, 1980) to UNIX time (January 1, 1970)
            if (nargout == 2)
                unix_time_f = rem(gps_time,1);
            end
            % constants in matlab are slower than copied values :-( switching to values
            % unix_time = uint32(gps_time - 0.5) + GPS_Time.UNIX_GPS_SEC_DIFF + uint32(gps_week) * GPS_Time.SEC_IN_WEEK;
            unix_time = uint32(gps_time - 0.5) + 315964800 + uint32(gps_week) * 604800;
        end
    end
    
    % =========================================================================
    %    GPS TIME
    % =========================================================================
    
    
    methods (Static, Access = 'public')
        % Shift from UNIX time (January 1, 1970 - msec) to GPS time (January 6, 1980 - sec)
        function [gps_week, gps_sow, gps_dow] = unixTimeToGps(unix_time, unix_time_f)
            % Conversion (shift) from Unix Time (January 1, 1970) to GPS time (January 6, 1980)
            
            % constants in matlab are slower than copied values :-( switching to values
            % gps_time = mod(double(unix_time -  GPS_Time.UNIX_GPS_SEC_DIFF, GPS_Time.SEC_IN_WEEK)) + unix_time_f;
            % gps_week = uint32(fix((unix_time - GPS_Time.UNIX_GPS_SEC_DIFF / GPS_Time.SEC_IN_WEEK));
            gps_sow = double(mod(unix_time - 315964800, 604800)) + unix_time_f;
            gps_week = uint32(fix(double(unix_time - 315964800) / 604800));
            gps_dow = uint32(fix(gps_sow / 86400));
        end
        
        function [gps_week, gps_sow, gps_dow] = date2gps(date_vec)
            % Conversion from calendar date to GPS time.
            % SYNTAX:
            %   [gps_week, gps_sow, gps_dow] = date2gps(date);
            %
            % INPUT:
            %   date = date [year, month, day, hour, min, sec]
            %
            % OUTPUT:
            %   gps_week = GPS week
            %   gps_sow  = GPS seconds of week
            %   gps_dow  = GPS day of week
            gps_start_datenum = 723186; %This is datenum([1980,1,6,0,0,0])
            
            %number of days since the beginning of GPS time
            %deltat   = (datenum([date(:,1), date(:,2), date(:,3)]) - gps_start_datenum);
            % hack: datenummmx is faster cause it does not check argins
            deltat   = (datenummx(date_vec(:,1:3)) - gps_start_datenum);
            
            gps_week = floor(deltat/7);            %GPS week
            gps_dow  = floor(deltat - gps_week*7); %GPS day of week
            gps_sow  = (deltat - gps_week*7)*86400;
            gps_sow = gps_sow + date_vec(:,4)*3600 + date_vec(:,5)*60 + date_vec(:,6); %GPS seconds of week
        end
        
    end
    
    % =========================================================================
    %    TESTS
    % =========================================================================
    
    methods (Static, Access = 'public')
        function test()
            % Testing function, tests some basic transformations
            log = Logger.getInstance(); % Handler to the log object
            
            log.addMessage('Testing Class Time - single value UTC');
            tic;
            t = GPS_Time(datenum('01-Jul-1994 00:00:00'), [], false);
            t.toGps();
            t.toUtc();
            t.toUnixTime();
            t.toMatlabTime()
            t.toRefTime();
            t.toMatlabTime()
            t.toUnixTime();
            t.toRefTime();
            t.toUnixTime();
            t.toMatlabTime();
            t_diff = abs(t.mat_time - datenum('01-Jul-1994 00:00:00')) * 86400;
            if t_diff < 1e-5
                log.addStatusOk('Passed: difference under threshold (0.01 ms)');
            else
                log.addWarning(sprintf('Difference greater than (0.01 ms): %e',t_diff));
            end
            toc
            
            log.addMessage('Testing Class Time - multiple consecutive values GPS');
            tic;
            mat_time = datenum('01-Jul-1994 00:00:00');
            mat_time = (mat_time(1)-1:1/86400:mat_time(1)+1)';
            t = GPS_Time(mat_time, [], true);
            t.toGps();
            t.toUtc();
            t.toUnixTime();
            t.toMatlabTime()
            t.toRefTime();
            t.toMatlabTime()
            t.toUnixTime();
            t.toRefTime();
            t.toUnixTime();
            t.toMatlabTime();
            t.toGps();
            t_diff = max(abs(t.mat_time - mat_time)) * 86400;
            if t_diff < 2e-5
                log.addStatusOk('Passed');
            else
                log.addWarning(sprintf('Difference greater than 0.2 ms: %e',t_diff));
            end
            toc
        end
    end
    
end
